Influence of the driving mechanism on the response of systems with athermal 
dynamics: the example of the random-field Ising model 
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We investigate the influence of the driving mechanism on the hysteretic response of systems with 
athermal dynamics. In the framework of local-mean field theory at finite temperature (but neglect- 
ing thermallly activated processes), we compare the rate-independent hysteresis loops obtained in 
the random field Ising model (RFIM) when controlling either the external magnetic field H or the 
extensive magnetization M. Two distinct behaviors are observed, depending on disorder strength. 
At large disorder, the //-driven and M-driven protocols yield identical hysteresis loops in the ther- 
modynamic limit. At low disorder, when the //-driven magnetization curve is discontinuous (due 
to the presence of a macroscopic avalanche), the M-driven loop is re-entrant while the induced field 
exhibits strong intermittent fluctuations and is only weakly self- averaging. The relevance of these 
results to the experimental observations in ferromagnetic materials, shape memory alloys, and other 
disordered systems is discussed. 

PACS numbers: 75.60.Ej, 75.50.Lk, 81.30.Kf, 81.40.Jj 



I. INTRODUCTION 

When slowly driven by an external force, systems with 
quenched- in disorder often respond via a sequence of ran- 
dom bursts or avalanches whose size distribution extends 
over several decades, ranging from microscopic to macro- 
scopic. A classic example is the Barkhausen noise that 
occurs in soft magnetic materials as one varies the ex- 
ternal magnetic field H. The random sequence of spikes 
observed in the induced voltage (that is, essentially, in the 
time-derivative of the induced magnetization M) is then 
related to the fluctuations in the motion of the magnetic 
domain walls due to structural disorder yj. Avalanche 
dynamics is also observed in other systems, as for in- 
stance in the acoustic emission associated to marten- 
sitic transformations or during the drainage of super- 
fluid 4 He from a mesoporous solid 3] . All these systems 
also exhibit rate-independent hysteresis when the exter- 
nal force is cycled up and down, thermal fluctuations be- 
ing negligible on the experimental time scale due to the 
presence of very high energy barriers (the systems then 
explore a set of long-lived metastable states) . 

Avalanches are characteristic of a system driven by an 
external force. In some experimental situations, however, 
it may be more convenient to control the conjugated ex- 
tensive variable (or its time-derivative) rather than the 
external field. In magnetic materials for instance, some 
suitable feedback control may be used so to impose the 
flux (or its time derivative) instead of the magnetic field. 
This allows to study materials that reverse their magne- 
tization very rapidly^! . In gas adsorption experiments, 
there are also cases where one controls the amount of 
gas introduced inside the porous material and measures 
the induced change in the pressure [5|- In shape mem- 
ory alloys, stress-strain curves are routinely obtained by 
controlling the deformation of the sample and measur- 



ing the stress [(j. The following question then naturally 
arises: what is the relation between the force-response 
diagrams (magnetization curves, adsorption isotherms, 
strain-stress curves, etc..) observed in the two situations, 
say "field-driven" and "magnetization-driven" (hereafter, 
we shall adopt the language of magnetic systems). This 
question is especially intriguing when there is a macro- 
scopic instability, e.g. a macroscopic jump in the mag- 
netization curve obtained with the field-driven protocol. 
In this case, one may force a system that would spon- 
taneously jump to a stable (or metastable) situation to 
remain into some other metastable or maybe unstable 
state. This may lead to a re-entrant behavior in the hys- 
teresis loop, as indeed observed experimentally jj, Q, 13 • 
To what extent this behavior reveals an intrinsic property 
of the free-energy landscape of the system or depends on 
the experimental protocol is not clear. To the best of our 
knowledge, this problem has not been seriously investi- 
gated in the literature, in spite of the intense theoretical 
activity focusing on out-of-equilibrium driven systems. 

In order to address this issue, we consider here a mi- 
croscopic model, the Random Field Ising Model (RFIM), 
which has been extensively studied at zero tempera- 
ture as a prototype of disordered hysteretic system with 
avalanche behavior In particular, in dimension 3, the 
RFIM is known to exhibit an out-of-equilibrium phase 
transition between a strong-disorder regime where the 
magnetization hysteresis loop is smooth on the macro- 
scopic scale and a low-disorder one where it displays a 
macroscopic jump. The nonequilibrium RFIM at T = 
is probably the simplest model to study the connection 
between metastability and hysteresis, and we expect its 
behavior in the magnetization-driven case to be relevant 
to a wide class of systems. However, as is discussed 
elsewhere [Toj. the definition of a metastable M-driven 
dynamics for Ising spins at T = is somewhat problem- 
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atic, and it is more convenient to work with continuous 
variables. This can be done in various ways and here 
we choose to study the RFIM at finite temperatures, us- 
ing a local mean-field approach in which each Ising spin 
is replaced by its thermal average (see e.g. Ref.[ll|). 
One then follows the evolution of the system among the 
local minima of the free-energy landscape, discarding 
all thermally activated processes that would allow for 
energy-barrier crossings. All the important features of 
the behavior of the RFIM at zero temperature (hystere- 
sis, avalanches, and out-of-equilibrium phase transition) 
are thus preserved. This approach is expected to be valid 
on an intermediate time scale during which the system 
can equilibrate locally but is not able to reach the global 
thermal equilibrium. In many disordered systems, the 
free-energy barriers are much larger than ksT and this 
intermediate time scale is the one probed by experiments. 

In the next section, we describe the model and briefly 
present the local mean-field approach in the field-driven 
case. The algorithmic method for the magnetization- 
driven situation is presented in section III and the nu- 
merical results in the strong and weak-disorder regimes 
are discussed in section IV. Some concluding remarks are 
given in section V. 



II. H-DRIVEN DYNAMICS FOR THE RFIM AT 
FINITE TEMPERATURE 

We consider the ferromagnetic (Gaussian) RFIM de- 
fined by the energy function 

u{{Si}) = -j SiSj - Yl hA (*) 

<ij> i 

where Si (i = 1, N) are Ising spins placed on the sites 
of a 3d cubic lattice of size Lx Lx L and hi are quenched 
random fields drawn from a Gaussian distribution with 
zero mean and standard deviation a. The first sum runs 
over all distinct pairs of nearest-neighbor sites. 

In the "unconstrained" situation that was considered 
in all previous studies^], the system is submitted to an 
external field H that is changed adiabatically. At zero 
temperature, it then tries to follow the field evolution by 
minimizing the magnetic enthalpy TL = Li — HM where 
M = Si is the overall magnetization. The relaxation 
dynamics generally used in this case is the T — version 
of the Glauber single-spin-flip dynamics where each spin 
is aligned to its local field /, + H, where 

/, •/»'• • ''• ( 2 ) 

j/i 

and the sum is over the nearest neighbors of site i (this 
corresponds to a local, partial minimization of 7i|l2j|'). A 
configuration {Si} is then (meta)stable according to the 
single-spin-flip dynamics when all the spins satisfy the 



stability condition 

5, = sign(/< + H) . (3) 

When H is changed from an initially stable state, a spin 
may become unstable, which may initiate an avalanche of 
other spin flips in the neighborhood. The avalanche ends 
when a new stable state is reached. The resulting hystere- 
sis loop (after averaging over disorder and extrapolating 
to the thermodynamic limit) is continuous for a > a c , 
with all avalanches being of microscopic size, and dis- 
continuous for a < a c , the discontinuity corresponding 
to a macroscopic avalanche that involves a finite fraction 
of the spins. The critical disorder for the RFIM on the 
cubic lattice has been numerically located at a c ~ 2.2 

urn. 

Extension of this behavior to finite temperature is 
straightforward when using a local mean-field approach. 
In this approximation, the Helmholtz free-energy of the 
system for a given realization of the random fields is a 
function of the thermally averaged local magnetizations 
mi, 

T{{mi\) = -J rniirij - h t mj 

<ij> i 

+ 28 + mi) ln( ^"^ } + (1 " mi) ln( ^"^ )} 

i 

(4) 

where (3 — 1/ (knT) |14|. In presence of an external field 
H, one considers the Gibbs free-energy 

g({m i },H)=F({m i })-HM (5) 

where M — X)i TO «- The minima of Q define the 
metastable states and satisfy the stationary conditions 
SQ/Suii — 0, i.e. 

rrn = tanh[/3(/i + H)} (6) 

where 

fi = J ^ m i + h% ( 7 ) 

j/i 

(in the following, the coupling constant J is taken as 
the energy unit and set to unity) . The method then con- 
sists in following the evolution of the system among these 
local minima as they are displaced or destroyed when 
the temperature T or the field H is changed [llj. Since 
one forbids any barrier crossing, the system is forced 
to stay in the same miminum until the stability limit 
is reached (which is signaled by the vanishing of some 
eigenvalue of the Hessian matrix d J- /dmidrrij). The 
system then jumps instantaneously to another nearby 
minimum. At constant temperature, this jump defines 
an iJ-driven avalanche. The dynamics has thus the same 
irreversible and hysteretic character as the T — single- 
spin-flip dvnamics|l5j|. The only effect of temperature is 
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to change the shape of the free-energy landscape. In par- 
ticular, there is only one minimum at high temperature 
and the evolution becomes reversible. 

In practice, Eqs. (6) and (7) are solved iteratively 



(i.e. m 



(n+l) 



tanh[/3(X; 



(n) 



hi + H)] where 



» 



is the value of mi after n iterations), changing H by 
small increments and using as initial configuration the 
m^'s obtained with the preceding value of H. A crucial 
property of such an iterative scheme is that it only con- 
verges to minima of the free-energy (and not to maxima 
or saddle-points) 0, . Thanks to the convexity of the 
hyperbolic tangent function in Eqs. (6), this "dynamics" 
has also the same properties as the T — dynamics: it 
is abelian[l8j (i.e. the metastable state reached after an 
avalanche does not depend upon the order in which the 
local magnetizations are updated during the avalanche) 
and it satisfies return-point memory |15|. 




FIG. 1: Ascending branch of the //-driven hysteresis loop 
computed within local mean-field theory at different tempera- 
tures in a system of size L = 20. The T = curve is obtained 
using the single-spin-flip algorithm as in Ref.[j]t. (Color on 
line) 

It is clear that the Gibbs free-energy obtained from Eq. 

(5) identifies with the exact enthalpy H = U — HM when 
T — > and rrii — > Si = ±1. The N coupled equations 

(6) then reduce to the stability conditions given by Eqs. 
(3). The iterative dynamics at finite temperature thus 
reduces to the usual single-spin-flip metastable dynam- 
ics when T — > 0. In particular, the hysteresis loop goes 
to the same T = loop, as illustrated in Fig. 1 below 
and above the critical disorder (in this figure and in the 
following, the temperature unit is J/ks- The mean-field 
critical temperature of the pure system is thus T c = 6). 
On the other hand, since thermally activated processes 
are completely neglected, the T = critical point where 



a macroscopic avalanche appears for the first time in 
the thermodynamic limit is now replaced by a full line 
a c (T) of out-of-equilibrium critical points [n|, as shown 
schematically in Fig. 2 (note that a precise numerical 
determination of the critical line would require an exten- 
sive finite-size scaling study as in Ref.pjj). Moreover, 
since the theory fully accounts for the disorder-induced 
fluctuations that play a dominant role in random-field 
systems, the critical behavior along the line is expected 
to be non-classical and identical to the one observed at 
T = Od. 
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FIG. 2: Schematic out-of-equilibrium phase diagram of the 
//-driven RFIM in the local mean-field approximation. The 
critical line cr c (T) separates the weak and strong disorder 
regimes characterized respectively by the presence or absence 
of a macroscopic jump in the hysteresis loop. (Color on line) 



III. M-DRIVEN DYNAMICS AT FINITE 
TEMPERATURE: ALGORITHM AND 
NUMERICAL METHOD 

We now consider the "constrained" situation where 
the magnetization of the system is changed externally 
and the field becomes an output variable. For simplic- 
ity, we assume that the magnetization m = M/N is 
monotonously increased from —1 to +1 (generalization 
to more complicated magnetization histories is straight- 
forward). In this case, the system tries to minimize its 
Helmholtz free-energy T instead of Q while satisfying the 
global constraint J2i m i = M . 

Since there is no intrinsic dynamics for Ising spins, 
there is no obvious way to define an M -driven dynamics 
for the RFIM that would be the counterpart of the T = 
single-spin-flip dynamics used in the //-driven case. In 
fact, it turns out that the simplest possible choice for a 
T = dynamics leads to an hysteretic behavior of the 
system that does not allow for dissipation, which does 
not seem very realistic^^- On the other hand, at finite 
temperature within the local mean-field approach, one is 
dealing with a standard problem of constrained optimiza- 
tion that can be solved in a natural way by the method 
of Lagrange mutipliers. One first introduces the function 
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C{{m i },\)=F[{m i })-\\^m i -M} (8) 

i 

where A is a Lagrange multiplier that has the meaning 
of a field coupled to the extensive variable J2i m i- Min- 
imizing J-({rrii}) with the constraint on the magnetiza- 
tion amounts to solve simultaneously the N + 1 coupled 
equations dC/drrii = and dC/dX = 0, 



m, = tanh^/i + A)] 



and 



M 



(9) 



(10) 



In a second step, one has to define an iterative scheme 
that tells the system how to go from one solution to the 
"nearest" one as the magnetization is slowly changed. 
For that purpose, we rewrite Eqs. (9) as 



(11) 



and sum over i so to express A as a function of m and 
the local magnetizations, 



A = — 1 1+m 

2/3 n ^I> 2 ^(l-rn,) 



(12) 



Eqs. (9), (7) and (12) define a set of nonlinear cou- 
pled equations that can be solved iteratively, as detailed 
below|20|. By construction, the configurations {to*}, so- 
lutions of Eqs. (9), are extrema of the Gibbs free-energy 
for the special value of the field A*(m) that satisfies the 
constraint equation, Eq.(10). In other words, the field 
A adjusts itself in the course of the iterations so that 
the system eventually reaches an extremum of Q with 
the correct magnetization. One can imagine that this 
is indeed the role of the feedback control used in mag- 
netic materials 0. It is worth noticing that the constraint 
y"V rrii = M is not satisfied until convergence is reached, 
which of course is the very principle of the Lagrange 
method. It is therefore unclear if one can attribute a 
physical meaning to the intermediate stages of the it- 
eration, although a similar problem may also occur in 
experiments [2l]]. 

In practice, to improve the convergence of the iteration 
procedure, we introduce a mixing parameter < a < 1 
that controls how much of the previous iteration is re- 
tained. When changing the magnetization from to to 
to + 8m, we thus search for a fixed point of the map 

m ( i n+1) = atanh[/3(jf > + A< n >)] + (1 - a)m l (n) 

A (n+1) = 0_ ^ 1 + TO 



+ (l-a)A (n) 
(13) 



-0,950 



-0,955 



ET -0,960 



-0,965 



-0,970 




•cc=0.95 
■cc=0.15 



0,4 0,6 0,i 



1,0 1,2 1,4 

(n) 



FIG. 3: Typical evolution of the system in the magnetization- 
field plane during the iteration procedure for two values of the 
mixing parameter a (cr = 1,T = 2, L = 20). m is varied from 
—0.97 to —0.96 (5m — 0.01). The inset is a magnification of 
the behavior during the last iterations. (Color on line) 



applied to all m, in parallel. Convergence is assumed 
when max(|m^" +1) - m\ n) \, 2/?|A(" +1 > - A^|) < e (with 
e = 1CP 4 in most of the calculations). 

The following general observations can be made: 
(i) In contrast with the -ff-driven situation, there is 
no obvious relation between the stability of the itera- 
tive scheme and the properties of the Hessian matrix |l7j. 
However, we strongly believe that convergence to a fixed 
point of the map (13) only occurs when the final con- 
figuration {to*} is a local minimum of Q, i.e. when 
d 2 J-/dmidmj is positive definite. This was checked in 
small systems (up to TV = 1000) by diagonalizing the 
Hessian for each solution {to* (to) , A* (to)}, using a stan- 
dard diagonalization routine 22] : in every case we found 
that all the eigenvalues were positive within numerical 
accuracy. In the following, we shall therefore assume 
that the configurations visited by the dynamics are al- 
ways metastable states. This, however, does not mean 
that these states are the same as those reached with the 
field driving, as will be discussed below. 

ii) Convergence does not occur at low temperatures 
(depending on a and L) . The iterative map either enters 
limit cycles that cannot be broken by further decreasing 
the value of a, or the constraint J2i m iW = M starts 
to be violated, whatever the value of e. This latter prob- 
lem is obviously related to the fact that some TOj be- 
come equal to ±1 within numerical accuracy. The prob- 
lem could perhaps be cured by working with transformed 
variables (see e.g. Ref.[23^ but we have not tested such 
a modified procedure. Of course, another possibility is 
that there are no more solutions to Eqs. (9), (7), and 
(12) at low temperatures. We shall come back to this 
problem in section V. 

iii) The Helmholtz free-energy J-({rrii}) may increase 
during the iteration procedure. This is in contrast with 
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FIG. 4: Portion of the H(m) curve computed in a sample 
of size L — 20 for a = 2 at T = 2, increasing m by very 
small increments 5m and changing the value of the mixing 
parameter a. The dashed (blue) curve is the corresponding 
m{H) curve obtained when driving the system with the field. 
(Color on line) 



what happens in the Jf-driven case where it is found 
empirically that G({mi}) always decreases during an 
avalanche (recall that H is kept constant during an 
avalanche). This is perhaps an indication that the tran- 
sient behavior described by the map (13) does not corre- 
spond to any reasonable physical dynamics on the free- 
energy landscape. But this may be a too pessimistic 
statement. 

iv) The field A varies in a complicated way during the 
iterations, as shown in Fig. 3. Accordingly, the evolu- 
tion of the local magnetizations with m is not necessarily 
monotonous and the numerical results may depend on the 
choice of the increment 5m and/or the mixing parameter 
a. In principle, this could be a serious shortcoming of 
the iterative scheme implying that the corresponding dy- 
namics is not deterministic (or "adiabatic" in the sense 
of Ref . [ill ) . i.e., the final state depends on how the map 
(12) is applied. Note for instance that the final values of 
A shown in the inset of Fig. 3 slightly differ for a = 0.95 
and a = 0.15. However, there is good numerical evidence 
that an adiabatic behavior is recovered when 5m is small 
enough. This is illustrated in Fig. 4 showing a portion of 
the curve H(m) computed in a single disorder realization 
of size L = 20 for two values of dm and a (in order to 
facilitate the comparison with the if -driven magnetiza- 
tion curve, we now use the same notation H for the field 
A* (m) and plot the curve H(m) with m on the vertical 
axis). One can see that the trajectory is independent of 
a for 8m = 10 -3 . Moreover, the same final states are 
obtained with 5m = 10~ 3 or 10~ 4 . Unfortunately, one 
cannot compute a whole curve in a reasonable time using 
such small values of 5m and most of the calculations pre- 
sented in the following were computed with 5m = 10~ 2 . 
With this value, the dynamics is not fully adiabatic but 
the variations when changing a are reasonably small, as 



illustrated in Fig. 3. 



IV. RESULTS 

The main results of our study are summarized in Fig. 5 
where we compare the m(H) and H(m) hysteresis loops 
obtained with the fZ-driven and Af-driven dynamics, re- 
spectively, for three different values of a at T = 2 (well 
below the mean-field critical temperature of the pure sys- 
tem). These calculations were performed in samples of 
linear size L = 60 or L = 80 [24|. One can readily no- 
tice two remarkable features: i) the Af-driven hystere- 
sis loops display strong fluctuations with to, especially 
at low a, and ii) the behavior is quite different in the 
strong and weak-disorder regimes characterized respec- 
tively by the absence or presence of a macroscopic (H- 
driven) avalanche. (Although we do not know the exact 
critical line cr c (T), the changes in the shape of the m(H) 
curve in Fig. 5 clearly indicate that 2 < a c (T = 2) < 3). 




FIG. 5: Comparison of the hysteresis loops obtained in a 
single sample with the _ff-driven (dashed line) and M-driven 
(solid line) algorithms for three values of a at T = 2. The 
size of the sample is L — 60 (for a = 3) or L — 80 (for a = 2 
and 1). (Color on line) 

Before discussing these results in more detail, let us 
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first comment on the mechanism that underlies the be- 
havior of the induced field as one slowly increases the 
magnetization. At the beginning, when to is close to —1, 
H has a large negative value and follows to smoothly 
because there exists only one stable state, the equilib- 
rium one. Hence the two curves H(m) and m(H) are 
identical. Then, for a certain value of m, there is a dis- 
continuous jump in the field towards a lower value. This 
corresponds to the first avalanche in the m{H) curve. 
One can thus associate this jump to the disappearance 
of the initial minimum: the system has to find another 
metastable state with the new magnetization. One then 
again follows this mimimum smoothly, until it reaches 
its own stability limit. Then, there is a new jump to the 
left and the evolution proceeds like this until to = +1. A 
similar behavior is observed along the descending branch. 
As can be seen in Fig. 4, when increasing to, the excur- 
sions of the field are always to the left of the iJ-driven 
magnetization curve and there may be several disconti- 
nuities in H{ni) corresponding to a single avalanche in 
m(H). Indeed, because of the ferromagnetic nature of 
the interactions, one can prove (at T = 0[25|], but also at 
finite temperature in the framework of local mean-field 
theory) that all the metastable states are inside the H- 
driven saturation loop. 

If this interpretation is correct, the magnitude of the 
horizontal jumps in H(m) is thus related to the orga- 
nization and the number of metastable states inside the 
-H-driven hysteresis loop. This explains the two different 
behaviors observed above and below the critical disorder 
<Jc{T). 



A. Strong-disorder regime 




FIG. 6: Evolution of the ascending branch of H(m) with sys- 
tem size for a — 3, 2 and 1 at T = 2. For the sake of clarity, 
the curves for the two larger sizes are displaced horizontally 
by 0.5 and 1. (Color on line) 



We first consider the strong-disorder regime a > a c (T). 
As can be seen in Fig. 5(a), the fluctuations in H{m) are 
small, which may be traced back to the fact that there are 
many metastable states in the vicinity of the iJ-driven 
magnetization curve. More precisely, one expects the 
number of metastable states inside the hysteresis loop 
to scale exponentially with TV when a > cr c (T). This 
was proven analytically in Ref.[2{| for the 1-d RFIM at 
T = and is probably always true in the strong-disorder 
regime where the saturation loop is continuous. In par- 
ticular, Fig. 4 in Ref.[25| shows that the complexity, 
i.e. the logarithm of the number of metastable states di- 
vided by N, increases very rapidly in the vicinity of the 
boundary (whereas it is zero outside the loop, indicating 
that there are no metastable states in this region when 
N — > oo). Since, moreover, the TJ-driven avalanches re- 
main of microscopic size when a > o~ c (T) , we expect the 
small wiggles in the H(m) curve to become infinitesi- 
mally small in the thermodynamic limit. This is indeed 
what is observed in Fig. 6(a) where the curves obtained 
for different system sizes are compared. Therefore, we 
predict that the two curves H(m) and m(H) should be- 
come identical in the thermodynamic limit, as suggested 



by Fig. 5(a). 

Note that these heuristic and numerical arguments im- 
ply that H(m) is a self-averaging quantity for a > a c (T). 
H{m) has a different value for each realization {hi} of the 
random fields and the issue of self-averaging concerns 
the behavior of the width of the probability distribu- 
tion Pl (if (to)) as L increases. If H(m) is self-averaging, 
"most" realizations (in a sense that can be made precise) 
lead to the same value of H(m) in the thermodynamic 
limit. If not, the measurement of the field in a single 
sample, no matter how large, does not give a meaningful 
result and must be repeated on many samples. The sit- 
uation is clearly different from the usual i/-driven case 
where the magnetization is the output variable: H(m) 
is indeed an intensive quantity and imposing the magne- 
tization M of a large sample does not impose the mag- 
netization in large subsamples. Therefore, the standard 
Brout argument 26J cannot be applied and one cannot 
say anything a priori about the variations of the variance 
Vff(m) = (H(m) 2 ) - (H{m)) 2 with L (hereafter, (...) de- 
notes the quenched average over disorder). However, if 
the two curves m(H) and H(m) are indeed the same in 
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FIG. 7: Variance Vfr(m) of H(m) for selected values of m as 
a function of system size for a = 3 and T — 2. The number of 
samples is 10000 for L — 20 and 1000 for L = 30, 40 and 60. 
For the sake of clarity, the variances for m = —0.96,-0.85 
and —0.75 are divided by 2, 4 and 6, respectively. The lines 
are fits to the form Vu{m) ~ L~ p , yielding p ~ 3 in all cases. 
(Color on line) 



the thermodynamic limit, H(m) is a self-averaging quan- 
tity, like m(H). This is confirmed by Fig. 7 that shows 
that the variance V#(m) does vary like L~ 3 (note that 
this calculation has only been performed in the initial 
part of the curve, in the range — 1 < m < —0.75, so to 
reduce the computational effort). 
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FIG. 8: Magnification of the lower part of Figs. 5(b) and 5(c). 
The dashed and dotted-dashed lines represent the metastable 
and unstable branches of the global mean-field curve obtained 
from Eq. (14). The arrows indicate the field H sp (a) at the 
stability limit. (Color on line) 



obtained by solving the implicit equation 27] 

r+ao 

(to) = / tanh\p(z{m)+h + H)]P(h)dh (14) 



B. Weak-disorder regime 

We now consider the case a < <r c (T) illustrated in Figs. 
5(b) and 5(c). There are two striking differences with Fig. 
5(a): i) the curve H(m) shows a pronounced re-entrant 
behavior, and ii) the field fluctuates much strongly with 
to. In line with the preceding discussion, we attribute 
the re-entrant behavior to the peculiar distribution of 
the metastable states in the field-magnetization plane in 
the weak-disorder regime. This issue was discussed in 
Ref.[25| where it was suggested (but not proven) that 
the macroscopic jump in the m(H) curve is related to 
the absence of metastable states in a certain region of the 
field-magnetization plane (see Fig. 10 in that reference). 
If this scenario is correct, the field, in the M-driven sit- 
uation, has to jump to the other side of this region so 
to find metastable states with the correct magnetization. 
We thus suspect that the H(m) curve is close to (or may 
be, in the thermodynamic limit, identical to) the bound- 
ary of this region. We want to stress that the re-entrant 
behavior discussed here has nothing to do with the S- 
shape predicted by the global mean-field theory where 
all the mj's are taken equal to the average magnetiza- 
tion (to). In Figs. 8(a) and 8(b), the lower part of the 
curves H(m) are compared to the mean-field prediction 



where z = 6 is the lattice coordination number and P(h) 
is the Gaussian probability distribution. As can be seen 
in Fig. 8(a) for er — 2, when increasing to, the field H(m) 
deviates from the mean-field curve much before the spin- 
odal H sp (a) is reached and it starts to strongly fluctuate 
with to just at the onset of the macroscopic avalanche 
in m(H) (note that the intermediate part of the mean- 
field curve is unstable whereas all the states visited by the 
A/-driven dynamics are metastable). Clearly, an average 
description cannot properly describe the metastable be- 
havior associated to the existence of local inhomogeneities 
in the system (it is only when to is close to ±1 and the 
system is quasi-homogeneous that the local and global 
descriptions become similar). 

The same is true for a = 1 in Figs. 5(c) or 8(b), 
but the re-entrant behavior is much more pronounced 
and the field almost drops to zero. Indeed, one ex- 
pects the metastable states to concentrate in the vicin- 
ity of H — as the disorder decreases (see for instance 
Fig. 1 in Ref. 25]). When a = 0, the only metastable 
states that remain in non-zero field are the two extremal 
ones that correspond to the metastable branches of the 
mean-field curve for the pure system, solution of m — 
tanh[/3(crn + H)]. Therefore, in the limit a — > 0, as to is 
increased from to = — 1, the field H(m) follows the lower 
metastable branch up to H sp (a = 0), then drops to zero 
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and stays to this value until the upper metastable branch 
is recovered and followed up to m = +1. The interme- 
diate metastable states in zero field then correspond to 
bubbles of the minority phase that grow progressively as 
m is increased (controlling m is like controlling the spin- 
odal decomposition of the system) . When a is finite but 
small, one expects a similar scenario: initially, the sys- 
tem is able to follow the imposed magnetization while 
remaining quasi-homogeneous; then, for a certain value 
of m (and a corresponding value of H that is smaller 
than H Bp (a) as can be seen in Figs. 8(a) and 8(b)), this 
situation becomes unstable and a nonequilibrium multi- 
domain state is formed while the field decreases suddenly; 
this multi-domain state then evolves progressively until 
the system becomes homogeneous again. This scenario 
may be compared with the one presented in Ref.0 (see 
p. 321 in that reference) to explain the re-entrant be- 
havior observed in picture-frame single crystals when the 
rate of change of the magnetic flux is imposed by using a 
feedback control psj. In this case, the motion of a single 
domain wall is the dominant process. However, one can 
notice the similarity of the curve plotted in Fig. 5(c) with 
the experimental curve psj. A similar behavior is also ob- 
served in strain-stress curves measured in shape memory 
alloys under displacement control loading condition (see 
e.g. Refs.0!)). In this case, the physical origin of the 
sudden drop in the stress after it has reached a critical 
value is still debated (in the recent simulations of Ref . [2^ , 
it is associated to the formation of a multi-domain state, 
like in the present work). Our results suggest that this 
is probably a generic behavior of the hysteresis loops ob- 
tained in weakly disordered systems when controlling the 
extensive variable. 

We now move to the discussion of the discontinuous 
fluctuations of the field H with the magnetization. As 
shown in Figs. 6(b) and 6(c), the magnitude of these fluc- 
tuations decrease with system size but they are still very 
significant for L — 100, the largest size investigated in the 
present work (we note however that the strongest discon- 
tinuities observed in Fig. 6(c) are displaced to larger and 
larger values of m as L increases, which suggests that 
they disappear as L — > oo). One may wonder whether 
one should observe such a sporadic, discontinuous behav- 
ior of the induced field in a macroscopic sample. This is 
again related to the issue of self-averaging. As shown 
in Figs. 9(a) and 9(b) (for selected values of m in the 
lower part of the H(m) curve), we find that the vari- 
ance Vff(m) in the weak-disorder regime still decreases 
with L, but the finite-size scaling exponent becomes sig- 
nificantly smaller than 3 (p ~ 1.75) as soon as m is in 
the range of the macroscopic avalanche. If this behavior 
does not change for larger L, this seems to indicate that 
H(m) is only weakly self-averaging j33, a rather remark- 
able behavior for a non-critical system. On the other 
hand, one may note that the corresponding histograms 
shown in Fig. 10 are still very well fitted by a Gaussian 
distribution. 

What happens just at the onset of the avalanche for 
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(a) o=2.0 



• m=-0.98 
■ m=-0.96 

♦ m=-0.85 
=-0.75 




FIG. 9: Same as Fig. 7 for a — 2 and a = 1. For the sake of 
clarity, the variances for a — 2 and m = —0.96, —0.85, —0.75 
are mutiplied by 2, 4 and 8, respectively. The variance for 
(j = l and m = —0.96 is not shown, as it does not decrease 
with L (see text) and the one for m = —0.75 is multiplied by 
3. The fit of the variance to the form Vjr(m) ~ L~ p yields p 
varying from 3 to ~ 1.75 as m increases from —0.98 to —0.75. 
(Color on line) 



(7 = 1 deserves a separate discussion. In this case, the 
sample to sample fluctuations of the field are so large 
that the variance does not decrease at all with the sys- 
tem size. Indeed, as can be seen in the right pannel of 
Fig. 10 for m = —0.96, H fluctuates between two quite 
distinct values, H max w 1.4 and H m i n * 0.4 (in Fig. 10, 
the two peaks are simultaneously present for L — 40 only, 
but a similar behavior occurs for other values of L and 
slightly different values of m). This double-peak struc- 
ture, which does not exist for a = 2, could indicate the 
occurence of a true discontinuity in H(m) in the thermo- 
dynamic limit, associated with the sudden appearance of 
the inhomogeneous multi-domain state. Rather, we be- 
lieve that the presence of the two peaks is a finite-size 
artifact: the number of metastable states between H m in 
and H max is probably so small that much larger statis- 
tics would be needed to observe intermediate values of 
the field. It is more likely that the drop in the field at 
the onset of the multi-domain state is very steep but that 
the curve H(m) is nevertheless continuous. 

A similar behavior is observed for the Helmholtz free- 
energy J 7 as a function of the magnetization: it is strongly 
self-averaging at large disorder (with a variance decaying 
like 1/i 3 ) and only weakly self- averaging when there is 
a macroscopic iJ-driven avalanche. In this latter case, 
however, the exponent p seems to be larger (p ~ 2.45) 
than the one associated to the fluctuations of the field. 
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physical behavior, as shown by the study of fluid adsorp- 
tion in porous solids 19]) but the dynamics is essentially 
a zero-temperature dynamics since all metastable states 
have an infinite life-time and free-energy barriers can- 
not be crossed. All the preceding calculations were thus 
performed rather arbitrarily at the reduced temperature 
T = 2 (below the mean-field critical temperature of the 
pure system), as we were mainly interested in the qual- 
itative behavior of the system as a function of disorder 
strength. The problem is that the limit T — > seems 
to be non-trivial. On the one hand, one would expect 
a gradual evolution with T, as is the case with the H- 
driven protocol (see Fig. 1). As an illustration, we show 
in Fig. 11 the evolution of the ascending branch of H(m) 
as T is decreased from 2 to 1. Not surprisingly, the fluc- 
tuations with m increase as T decreases but one could 
imagine that it is possible to decrease the temperature 
down to zero. This optimistic view is reinforced by the 
hope that the envelope of the metastable states in the 
H ~ M plane (which we think plays an important role in 
understanding the behavior of the M-driven system) has 
a well-defined limit at T = 0, as discussed in Ref . |25| . 



FIG. 10: Histograms of H(m) for different values of m in the 
range — 1 < m < —0.75 for a = 1,T = 2 and different system 
sizes. The statistics is performed over 10000, 1000, 1000 and 
600 disorder realizations for L = 20, 30, 40 and 60, respec- 
tively. (Color on line) 



V. CONCLUDING REMARKS 

To summarize, we have presented a detailed study of 
the nonequilibrium response of the random field Ising 
model when one cycles adiabatically the magnetization 
(M-driven protocol) instead of the magnetic field (H- 
driven protocol). The study has been performed in the 
framework of the local mean-field approach at finite tem- 
perature, neglecting all thermally activated processes. 
We have shown that two regimes can be observed, de- 
pending on the disorder strength. In the strong disorder 
regime, where all .ff-driven avalanches are of microscopic 
size, the two protocols yield the same hysteresis loop in 
the thermodynamic limit. On the other hand, in the 
weak-disorder regime, where the iJ-driven magnetization 
curve is discontinuous, the M-driven hysteresis loop is 
re-entrant and the output field displays strong fluctua- 
tions with the magnetization. There are also considerable 
sample-to-sample fluctuations that decrease only weakly 
with system size. 

From the results that have been presented above, a 
question naturally arises: what is the behavior of the 
system as the temperature is decreased to zero ? As 
noted in the introduction, temperature was introduced in 
the model so to work with continuous variables instead 
of Ising spins and to use the Lagrange constrained op- 
timization method. The free-energy landscape changes 
with T (which may be very important for the actual 
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FIG. 11: Effect of temperature of the ascending branch of 
H(m). The system size is L — 40. (Color on line) 

Unfortunately, the iteration procedure does not con- 
verge any more at lower T for reasons that are not com- 
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pletely clear. As noted earlier, this may be a problem of 
numerical accuracy in the solution of the iterative map. 
This may also be related to the fact that the magneti- 
zation M becomes a discrete variable when T — ► and 
can only change by an increment of ±2. It may be then 
very difficult for the algorithm to find a path from one 
metastable state to another due to the existence of these 
gaps in the magnetization. Clearly, it would be much 
more convenient to define a proper dynamics directly at 
T = 0. However, as discussed elsewhere [Toj . this does 
not seem to be an easy task either. 

Finally, we would like to remark that some of the fea- 
tures found in the present calculations, in particular the 
re-entrant behavior and the intermittent fluctuations in 
the induced field as a function of magnetization, are ob- 
served in a variety of systems characterized by an ather- 
mal dynamics. For instance, many amorphous materials, 
including polymers and metallic alloys, exhibit a mechan- 
ical response to deformation (e.g. the stress induced by 
a quasistatic shear strain) that closely resemble the field- 
magnetization curves computed in this work. Similar be- 
havior is also observed in granular materials and foams. 



The relation of this behavior to changes in the underly- 
ing energy landscape has been discussed in the literature 
(see e.g. Ref|3l|l. However, it seems to be the first time 
that a direct comparison between the different driving 
mechanisms is performed. 
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